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Abstract 

Gravitational and hydrodynamical perturbations are analysed in a relativistic plasma con- 
taining a mixture of interacting fluids characterized by a non-negligible bulk viscosity co- 
efficient. The energy-momentum transfer between the cosmological fluids, as well as the 
fluctuations of the bulk viscosity coefficients, are analyzed simultaneously with the aim of 
deriving a generalized set of evolution equations for the entropy and curvature fluctuations. 
For typical length scales larger than the Hubble radius, the fluctuations of the bulk viscosity 
coefficients and of the decay rate provide source terms for the evolution of both the curvature 
and the entropy fluctuations. According to the functional dependence of the bulk viscosity 
coefficient on the energy densities of the fluids composing the system, the mixing of entropy 
and curvature perturbations is scrutinized both analytically and numerically. 
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If a relativistic plasma contains a mixture of inviscid fluids with negligible transfer of 
energy and momentum, the evolution of entropy fluctuations is characterized, in the long- 
wavelength limit, by the absence of source terms containing curvature perturbations. This 
property has relevant consequences on the dynamics of the coupled system of gravitational 
and hydrodynamical perturbations. It implies, for instance, that curvature perturbations are 
conserved, in the long- wavelength limit, under rather general assumptions [1, 2, 3]. Long- 
wavelength fluctuations in the spatial curvature determine, via the Sachs- Wolfe effect, the 
temperature inhomogeneities observed in the microwave sky (see, for instance, [4]). 

Mixtures of relativistic fluids are a useful toy model that can be investigated with the pur- 
pose of inferring some general properties of the evolution equations of curvature and entropy 
fluctuations. Moreover, multifluid systems are per se relevant to the model-independent 
discussion of initial conditions of CMB anisotropy [1, 4]: for instance, the five isocurva- 
ture modes supported by the predecoupling plasma may be discussed, in their simplest 
realization, by a truncated Einstein-Boltzmann system of equations whose lower multipole 
moments reproduce indeed a multifluid hydrodynamical description [5]. 

One of the assumptions often invoked in the analysis of multifluid systems is that the 
bulk viscosity coefficient and its possible spatial variation have a negligible impact on the 
dynamics. While this assumption may be justified in some specific system, it may not be 
true in the early stages of the life of the Universe (see, for instance, [6, 7, 8]). Unlike other 
dissipative effects, the presence of bulk viscosity does not spoil the isotropy of the background 
geometry. Therefore, consider a mixture of two relativistic fluids (the a-fluid and the b-fluid) 
obeying a set of generally covariant evolution equations formed by the Einstein equations 2 

r; - = \v (!) 

and by the evolution equations of the energy-momentum tensors of each fluid of the mixture, 
i.e. 

v M rr = -r«f%( Pa + Pa ), (2) 
v M Tr = r^%( Pa + Pa ), (3) 

where up is the total velocity field of the mixture. Equations (2) and (3) describe the 
situation where the a-fluid decays into the b-fluid with decay rate Y. It is evident from 
the form of Eqs. (2) and (3) that the total energy-momentum tensor of the mixture, i.e. 
T^ v = T£ v + T£ v is covariantly conserved, i.e. V M T^ = 0. The total energy-momentum 
tensor of each species is given by the sum of an inviscid contribution, denoted by T^ h and 

2 Units of 87rG = 1 will be used throughout. Notice, to avoid confusions, that the Latin (lower-case 
roman) subscripts a, b, c d, ... will denote, in the present paper, different fluids present in the relativistic 
plasma. Greek (lower-case) subscripts will denote tensor indices. Latin (lower-case italic) subscripts i, j, k, ... 
will denote the spatial components of a tensor. 
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by a viscous contribution, denoted by T^ b , i.e. 

= ^b + Tfb, (4) 
T Ch = (Pa, b + Pa, b)< b «a, b ~ Pa, b^, (5) 

Tb = ea,b(^-< b < >b )V Q <b, (6) 

where the subscript in the various fluid quantities simply means that Eqs. (4), (5) and (6) 
hold, independently, for the a- and b-fluids. So, for instance, in Eqs. (5) and (6), u£ and ii b 
denote the peculiar velocities of each fluid of the mixture. 

In a spatially flat metric of Friedmann-Robertson- Walker (FRW) type charaterized by 
a background line element 



ds 2 = g^dx^dx" = a 2 (r)[dT 2 - dx\ (7) 



Eqs. (2) and (3) imply 



p' a + m{ Pa + v & ) + ar( Pa + Pa ) = o, (8) 

p' h + 3H(p b + Pb) - af(p a + p a ) = 0, (9) 

where the prime denotes a derivation with respect to the conformal time coordinate r and 
7i = a' I a. In Eqs. (8) and (9), 7\b denote the total effective pressure of each species, i.e. 

Pa,b=Pa,b-3^a,b, (10) 

while p a b denote the inviscid pressures of each species. In Eq. (10), £ b denote the bulk 
viscosity coefficient evaluated on the background geometry. As will be discussed later, the 
bulk viscosity coefficient may depend on both p a and Equations (8) and (9) lead to the 
evolution of the total energy and pressure densities 

p' + 3H(p + V) = 0, (11) 

where p = p a + Pb and V = V a + Vb- Equations (8) and (9) must be supplemented by the 
explicit background form of Eq. (1), i.e. 

m 2 = a 2 p, (12) 
2{H 2 -H') = a 2 {p + V), (13) 

where, again, (p + V) is the total effective enthalpy that contains the background viscosity 
coefficient of the mixture £ = £ a + £ b . 

We are now interested in deriving the evolution of the entropy and total-curvature fluc- 
tuations of the system. Both the entropy perturbations and the perturbations in the total 
spatial curvature can be written in terms of £ a and which are related in the off-diagonal 
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gauge [9] (see also [10, 11]) to the density contrasts of the individual fluids of the mixture 
[12]: 



S = -3(C.-Cb), ( 14 ) 

C = 4Ca+4Cb- (15) 

p p 

In the following we are going to exploit the off-diagonal gauge [9] (also called uniform- 
curvature gauge [10]), which is particularly convenient for the problem at hand. The results 
will be exactly the same of those obtainable in the framework of gauge-independent descrip- 
tions (see [13]). In fact, the quantities S and (, defined in terms of £ a and are invariant 
under infinitesimal coordinate transformations. Consequently they can be computed in any 
suitable (non-singular) coordinate system. 

Concerning Eq. (15) it could be noticed that the gauge- invariant definition of the spatial 
curvature perturbation can be slightly different for wavelengths smaller than the Hubble 
radius. However, since in the problem at hand we are mainly interested in super-Hubble 
fluctuations, Eq. (15) is numerically equivalent to the curvature inhomogeneities defined, 
more conventionally, from the curvature fluctuations on comoving orthogonal hypersurfaces 
[2, 14, 15] (see also [4]). The physical interpretation of the entropy fluctuations defined in 
Eq. (14) can be understood, for instance, in the case of an inviscid mixture of cold dark 
matter (CDM) and radiation. In this situation the a-fluid is given by pressureless matter 
and the b-fluid by radiation. Applying the covariant conservation equations for two (inviscid 
and non-interacting) species, it is clear that S is nothing but the fractional fluctuation in 
the specific entropy q = T 3 /n CDM (where T is the radiation temperature and ucdm is the 
CDM number density), i.e. the entropy density per CDM particle: 

S = — = Scdm - 7<5 r , (16) 
q 4 

where <5cdm and S T are the density contrasts for the CDM and for the radiation fluids. 

In the off-diagonal gauge the spatial components of the perturbed metric vanish and, 
hence, the only components of the perturbed line element are: 

5g Q0 = 2a 2 0, Sgoi = -a 2 d { B. (17) 

Since, in the long-wavelength limit, (H' — H 2 )(f) = H 2 (, it turns out that in the off-diagonal 
gauge, 5goo is connected to (. As anticipated, ( a and can be expressed in terms of the 
fluctuations of the density contrasts of the individual fluids, i.e. 

Ca = ^a, Cb = —Ob- (18) 

Pa Pb 

The evolution equations obeyed by the density contrasts 5 a = Sp a /p a and 5^ = 5pb/Pb are 
derived by perturbing Eqs. (2) and (3) to first order in the amplitude of the metric and 
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hydrodynamical fluctuations: 



9H 2 - 

5' a + (3H + aT)(c 2 sa - w a )b a + [£ a (0 + Q - + a(l + w a )T(b r + 0) 



+ 



1 + 0-6 



e* = o, 

;Pa 



b' h + 3H(c 2 sh - w h )5 h + ar^[(l + w a )(b h - S T - 0) - (1 + cfJ5 a 

Pb 



ap h 



;i+^b)-6 



#h = 0. 



(19) 



(20) 



Concerning Eqs. (19) and (20) a few comments are in order: 

• for notational convenience the barotropic indices (i.e. w a , w^) and the sound speeds 
(i.e. c 2 a and c 2 b ) have been introduced for the inviscid component of each species of 
the plasma; if the inviscid component is parametrized in terms of a perfect relativistic 
fluid 4,b EW a,b; 

• b r = bT/T is the fractional fluctuation of the decay rate computed in the off-diagonal 
gauge; 

• 5£ a and 5£b denote the fluctuations of the bulk viscosity coefficients; later on it will also 
be convenient to introduce the fluctuation in the total viscosity, i.e. <5£ = <5£ a + b^; 

• finally, 9 a = <9jt> a = V 2 t> a and 8^ = <9jt> b = V 2 t>b are the divergences of the peculiar 
velocity field of each species; note that the global velocity 9 = diV 1 field (with 8Tq = 
(p+V)v l ) is recovered from 9 a and 9 h by recalling that (p+p)9 = (p a +Pa)# a +(Pb+Pb)#b- 

Equations (19) and (20) must be supplemented by the perturbed components of Eq. (1); 
in particular by the Hamiltonian and momentum constraints: 



HV 2 B + 3ft 2 + — bp = 0, 
2 

a 2 

V 2 [H(f> + (H 2 - H')B\ + — (p + V)9 = 0, 



(21) 
(22) 



and by the other two equations stemming from the spatial components (i.e. [i = j) and 
(i^))ofEq. (1): 



0' + (ft + 2^-- 
B' + 2KB + = 0. 



5 P -3— -m<f) 

a a 



0. 



(23) 
(24) 



In Eq. (21) the global energy and pressure density fluctuations (i.e. bp and bp) have been 
introduced. As is clear from Eqs. (21)-(24), one of the advantages of the off-diagonal 
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formulation is the absence of second time derivatives of the metric fluctuations. Strictly 
speaking the evolution equations of 6 a and 9^ should be added to the system. However, they 
are only relevant for typical length scales smaller than the Hubble radius at a given time. 
Since we are interested in the opposite regime, these equations will be omitted, but they will 
be discussed elsewhere in their full generality (see [13]). 

By combining Eqs. (21) and (23) the evolution equation for ( can be easily obtained; it 
is given by 



C = ~^[(L + £b)C + # (*& + %)] + fHc4 - c4)(Ca - C), 



2H 



2H 



(25) 



where we passed, for later convenience, from the conformal time coordinate r to the cosmic 
time coordinate t (i.e. dt = a(r)dr). 

Equations (19) and (20) lead to the evolution equations of Ca and Cb whose explicit form 
is given by 



Ca + 

and by 



+ (3# + r)(l + c; 



say 



5a 



Pa + Pa p 
5a 



6r + §- 2 C 



(26) 



Cb + 



% 



T?a /-, , 2 \ /■ , ^ r i 



C b - r^(i + C s 2 JCa + -[H^ h - H^Q 

5b 5b 



Pa + Pa p 
5b 



H , 

Sr + 



where 



Pa Pb 
5a = 77, 5b 



(27) 



(28) 



H 1 H 

Various identities can be used to bring Eqs. (25), (26) and (27) to slightly different (but 
equivalent) forms. In particular: 

• using Eq. (15), we can always trade the combinations (Ca~ C) an d (Cb — C) f° r Pa/p(Ca — 
Cb) and p b /p(Ca-Cb); 

• according to Eq. (14), (£ a - Cb) = S/3; 

• by virtue of the background equations (12) and (13), H/H = p/(2p); 

• if the inviscid component of each fluid of the mixture is a perfect fluid, then c 2 a = w a 
and c 2 b = w h ; 

• finally the background evolution of each fluid, i.e. Eqs. (8) and (9), may always be 
employed to obtain equivalent forms of the above equations. 
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Specific limits of Eqs. (25)-(27) will now be reproduced. In the limit £ a = £ b = 0, with 
§ r = 5^ = 5£ h = and f = 0, Eqs. (25)-(27) read 




In the case u> a = and = 1/3, Eqs. (29)-(31) coincide with the set of equations used in 
Ref. [16] to describe the radiative decay of a massive curvaton whose effective pressure, at 
the oscillatory stage, reproduces that of dusty matter, i.e. w a = 0. It is then clear, taking the 
difference between Eqs. (31) and (30), that the evolution equation of entropy perturbations 



is homogeneous and does not contain any ^-dependent source term. 

Sticking to the case of the radiative decay of a dusty fluid, but including the fluctuations 
of the decay rate, the following system of evolution equations 



can be derived from Eqs. (25)-(27) when w\> = 1/3 and u> a = 0. In Eqs. (33) and (34), 
g a = — Hp a /p a . Equations (33) and (34) describe the situation discussed in Ref. [17], where 
the dynamics of the inflaton with inhomogeneous decay rate has been discussed (see, for 
instance, also [18, 19, 21, 20]). If the spatial fluctuations of the decay rate are not a function 
of the local energy density of the mixture, curvature fluctuations may be generated for length 
scales larger than the Hubble radius. 

Consider now the case where the T is constant, the decay is homogeneous (i.e. 5? — 0), 
but £ a = Ca(Pa) and £b = 6t>(Pb)- This occurrence implies that 




(32) 



C = ^(Ca-C), 



(33) 




(34) 




(35) 




■) 



S, 



(36) 



(37) 



(38) 
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Again, in this case, it can be easily argued that the evolution of entropy fluctuations obeys 
a homogeneous equation in S. In fact, combining Eqs. (37) and (38) it is possible to obtain: 



S = - 



!^rK + i)(i-4-2^ 

2pPb V pi p a 



"a 



s. 



(39) 



This conclusion can be, however, evaded if £ a and £b are functions both of p a and p b , i-e. 
£a = £ a (p a ,Pb) and £ b = Cb(p a ,Pb)- In this case 



= -|(Ca + Cb), 56 = -|(Ca + Cb)- 



H 



(40) 



Thus, in the situation described by Eq. (40), Eq. (39) will inherit two extra terms at the 
right-hand side, i.e. 

9 H 2 ( - - \ 

— PaCbCa - PbCaCb , (41) 
Pb Pa V / 

which cannot be recast, for generic £ a and £ b , in a single term proportional to S. 

A relevant issue to be addressed concerns the phenomenological viability of interacting 
viscous mixtures. Consider, for instance, a model where the decay rate is constant but 
inhomogeneous (i.e. 5r ^ 0) and £ a = a^fpl (where e is constant). The viscosity coefficient 
of the b-fluid vanishes, i.e. £ b — 0. This model describes the situation where the a-fluid is 
initially dominant and characterized by a viscosity proportional to e. Furthermore, if we want 
the Universe to be expanding, we must also require e < (w a + l)/y/3. The a-fluid will start 
its decay for a typical cosmic time tr ~ T , and then the background will be dominated 
by the b-fluid while the energy density of the a-fluid, i.e. p a , will decay exponentially. 
Also the background viscosity will decay exponentially, since £ a = e^/p^. These aspects are 
illustrated in Fig. 1 (plot at the left-hand side) where, for two different values of e, the 
common logarithm (i.e. the logarithm to base 10) of p a and pb are reported. From the point 
of view of the background, this model is perfectly viable and it leads to a final stage of 
expansion dominated by the b-fluid. To make the example even more explicit, one can think 
of the situation where the a-fluid is given by dust (i.e. u> a = 0) or stiff matter (i.e. w a = 1). 
The b-fluid may be taken, for instance, to coincide with radiation (i.e. u>b = 1/3). 

The dynamics of curvature fluctuations may be described, for practical reasons, by ex- 
pressing the evolution equations in terms of ( and 5 a , i.e. the curvature fluctuations and the 
density contrast of the a-fluid. Given the relations (15) and (18) Cb 5 Ca and 5 b can always be 
obtained as linear combinations (with background-dependent coefficients) of ( and 5 a . From 
Eqs. (19) and (25) the relevant evolution equations can be written as 



c = 



AH 



3^— + 2(w h - w a 



(paC + H p a S a ), 



• 9eH\ H 
5 a + --=5 a + 



2^ 



H 2 



9eH 2 

VPa" 



+ I> a + 1) 



C = -i> a + i)<$r. 



(42) 
(43) 
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Figure 1: The evolution of the background (left-hand plot) and of the fluctuations (right- 
hand plot) is illustrated. The parameters of the mixture are fixed in such a way that w a = 0, 
tfb = 1/3, T/Hq ~ 10~ 3 , £ a = e^/p^, £b = 0; H denotes the value of the Hubble parameter at 
the initial integration time. In the left plot, the dashed curve represents the evolution of the 
energy density of the decay products (radiation) while the full lines represent the evolution 
of the energy density of the decaying component for different values of e. In the right plot 
hand side the dashed curves illustrate the behaviour of |£ a | while the full lines represent the 
evolution of \(\. Both \(\ and |£ a | are given in units of Sr- 



Equations (42) and (43) describe the evolution of ( and <5 a for typical wavelengths larger than 
the Hubble radius. Initial conditions of the system are then set by requiring C(^o) — and 
^a(^o) — $b — 0) where to is the initial integration time. From Fig. 1 (plot at the right-hand 
side) the evolution is such that curvature fluctuations grow from to a value proportional 
to Sr, i.e. proportional to the fluctuations of the decay rate over length scales larger than 
the Hubble radius. The final asymptotic value of ( can be determined analytically and it 
turns out to be 

1 (\ + 3 v / 3e^ 



In the limit e — > 0, the results reproduce the findings of Ref. [17] leading to a Bardeen 
potential |^fi na i| — <5r/9, which implies |Cfinai| — by using the well-known relation of ( 
and \1/ in a radiation-dominated phase (see, for instance, [4]). In the example discussed so 
far the values of S-p have been taken in the ranges 10~ 6 -10~ 9 . 

The class of examples reported so far can be generalized in various ways. Different 
barotropic indices for the fluids of the mixture can be studied. Equation (44) can then be 
generalized to the cases of generic w & and w^. Furthermore, the functional dependence of the 
viscosity coefficients can be chosen to be different. Possible generalizations will be present 
elsewhere [13]. We would like to point out that the simple examples presented here may be 
made more realistic by thinking that a dust fluid is an effective description of a scalar field 
oscillating in a quadratic potential [22]. Thus, the simple fluid model of a dust fluid decaying 



into radiation has been used [17] (with some caveats [20]) to infer some properties of the 
infiaton decay when the inflaton decay rate is not homogeneous. If the inhomogeneous decay 
occurs after an inflationary phase at low curvature (i.e. H in { <C lCT 6 Mp), it is plausible to 
argue that the spectrum of 5r may be converted into the spectrum of ( for typical frequencies 
smaller than the Hubble rate. We are not interested here in supporting a specific model of 
inhomogeneous reheating. The purpose of the examples discussed so far is purely illustrative. 
However, the lesson to be drawn is that bulk viscous stresses may play a relevant role. 

In the present paper, various results have been achieved. First of all, the concept of 
interacting viscous mixtures has been introduced, i.e. a mixture of interacting fluids with 
viscous corrections. In this framework, the coupled evolution of curvature and entropy 
fluctuations has been derived in the case where both the decay rate and the bulk viscosity 
coefficients are allowed to fluctuate over typical length scales larger than the Hubble radius. 
Different situations have been systematically discussed. If the decay rate is constant and 
homogeneous, with bulk viscosities that depend separately on the energy density of each fluid 
of the mixture, the evolution of entropy fluctuations obeys a source-free evolution equation. 
If, on the contrary, the bulk viscosity has a more general dependence on the energy densities 
of the fluids composing the mixture, the evolution equations of the entropy perturbations 
may inherit a source term that involves, in one way or another, curvature fluctuations. In 
similar terms, if the decay rate is allowed to fluctuate without being a function of the local 
density of the fluid, entropy fluctuations will not obey a source-free equation. 
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